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ABSTRACT 


Random  Wave  Forces  on  a  Free-to-Surge  Vertical  Cylinder  (May  1988) 
Charles  Blake  Sajonia,  B.S.,  University  of  Washington 
Chairman  of  Advisory  Committee:  Dr.  J.  M.  Niedzwecki 


The  principal  objective  of  this  research  is  to  to  gain  insight  into  the  applications 
and  limitations  of  the  relative  motion  form  of  the  Morison  equation  for  the 
prediction  of  hydrodynamic  forces  on  a  free-to-surge  vertical  cylinder  in  random 
waves.  Force  transfer  coefficients  are  estimated  from  experimental  data  using 
regression  with  autoregressive  errors.  The  best-fit  relative  motion  form  of  the 
Morison  equation  results  in  a  root-mean-square  error  of  24%  and  a  multiple 
correlation  coefficient  of  0.71,  over  a  16  second  time  series.  A  high  frequency  force 
component  not  accounted  for  in  the  Morison  equation  is  quantified.  Cross-spectra 
are  used  to  show  that  this  residual  force  can  not  be  duplicated  by  the  relative  motion 
Morison  equation  due  to  the  lack  of  explicit  history  terms.  A  force  transfer  model 
containing  explicit  history  terms  is  presented.  The  improvement  in  force  prediction 
with  increasing  memory  is  illustrated  and  a  memory  length  is  chosen  that  optimizes 
the  tradeoff  between  model  complexity  and  goodness-of-fit.  The  new  model  reduces 
the  rms  error  from  24%  to  9%,  increases  the  multiple  correlation  coefficient  from 
0.71  to  0.83,  and  captures  the  high  frequency  force  components  not  accounted  for 
in  the  Morison  equation.  A  simple  numerical  simulation  of  a  tension  leg  platform 
is  performed  to  illustrate  the  application  and  limitations  of  the  results. 
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CHAPTER  I 
INTRODUCTION 


1 


A.  The  Relative  Motion  Morison  Equation 


The  principle  objective  of  this  research  is  to  gain  insight  into  the  applications 
and  limitations  of  the  relative  motion  form  of  the  Morison  equation  for  the 
prediction  of  hydrodynamic  forces  on  a  free-to-surge  vertical  cylinder  in  random 
waves.  In  1950  Morison,  O’Brien,  Johnson,  and  Schaaf  (1950)  presented  an 
empirical  equation  which  describes  the  wave  induced  hydrodynamic  loading  on 
a  fixed  vertical  pile.  Their  equation,  which  is  popularly  known  as  the  Morison 
equation,  can  be  expressed  as 


,  _  D  ,  ,  „  irD2 

J  =  Cdp—\u\u  +  Cmp  —  a, 


(1.1) 


where  /  is  the  force  per  unit  length  acting  on  the  pile,  Cd.  and  Cm  are  force  transfer 
coefficients,  u  and  a  are  the  horizontal  fluid  particle  velocity  and  acceleration,  p  is 
the  density  of  the  fluid,  and  D  is  the  diameter  of  the  pile  or  cylinder.  The  total 
wave  induced  force  acting  on  the  cylinder  is  considered  as  the  sum  of  a  viscous  drag 
force  component  and  an  inertia  component.  When  the  structure  moves  in  response 
to  hydrodynamic  loads,  Equation  (1.1)  is  often  modified  by  using  relative  velocities 
and  accelerations.  The  resulting  equation  is  known  as  the  relative  motion  or  wave- 
structure  interaction  form  of  the  Morison  equation.  This  equation  can  be  expressed 
as 


,  _  D  ...  x  „  ttD2, 

/  =  Cdp  —  \u  -  x|(u  -  i)  +  Cap-j-{a 


7 vD2 

x)  +  P-~r~a > 


(1.2) 


Journal  model  is  ASCE  Proceedings  of  the  Waterway,  Port,  Coastal  and  Ocean 
Engineering  Division. 
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where  x  and  x  are  the  velocity  and  acceleration  of  the  cylinder,  respectively,  Ca  is 
the  added  mass  coefficient  defined  as  Ca  =  Cm  ~  1?  and  the  last  term  is  known  as 
the  Froude-Krylov  force  which  is  the  result  of  the  local  acceleration  of  the  unsteady 
flow. 

A  number  of  researchers,  e.g.  (Sarpkava  and  Isaacson,  1981),  have  shown  that 
the  force  transfer  coefficients  Cd  and  Cm  are  not  simple  constants,  but  functions  of  at 
least  four  parameters:  cylinder  roughness,  Reynolds  number,  Keulegan-Carpenter 
number,  and  time.  Thus,  the  relative  motion  Morison  equation  is  an  engineering 
approximation  to  a  complex  problem. 

B.  Flow  History 

Although  the  Morison  equation  contains  only  instantaneous  or  time  averaged 
parameters,  the  actual  hydrodynamic  force  acting  on  a  submerged  body  depends  on 
both  the  instantaneous  and  preceding  flow  conditions  (Hamiliton,  1972).  However, 
most  investigators  have  ignored  this  history  dependence  in  light  of  the  many 
other  uncertainties  in  wave  force  calculations.  Others,  such  as  Keulegan  and 
Carpenter  (1958)  and  Sarpkaya  (1981)  have  corrected  for  this  error  by  including  an 
instantaneous  correction  term.  In  periodic  flow,  the  fluid  particle  kinematics  change 
with  time,  but  the  pattern  is  continually  repeated.  Since  the  preceding  cycles  are 
identical,  the  history  effects  can  be  represented  by  an  instantaneous  correction  term 
instead  of  an  explicit  history  term.  This  approach  was  utilized  by  Keulegan  and 
Carpenter  (1958)  in  a  series  ol  seiche-tunk  flow  experiments.  A  “remainder  force” 
A R  was  introduced  to  account  for  the  differences  between  predicted  and  measured 
forces  on  a  submerged  cylinder. 

Sarpkaya  (1981)  has  also  quantified  A R  in  a  series  of  one-dimensional  periodic 


flow  experiments.  In  addition,  the  limitations  of  the  Morison  equation  were 
illustrated  by  plotting  the  instantaneous  values  of  Cd  and  Cm.  It  was  shown  that 
Cd  and  6  n  exhibit  large  variations  during  a  given  cycle,  particularly  for  values 
of  the  Keulegan-Carpenter  number  between  8  and  25,  and  that  accelerating  and 
decelerating  flows,  with  identical  absolute  values  of  the  corresponding  fluid  particle 
kinematics,  do  not  exert  identical  forces  on  the  cylinder.  These  observations  suggest 
that  A R  is  partly  the  result  of  neglecting  the  effects  of  flow  history  and  partly  due 
to  the  other  simplifications  inherent  in  the  Morison  equation.  Sarpkaya  represented 
the  ‘‘remainder  force”  as  an  instantaneous  correction  term  given  by  the  Fourier 
expansion  proposed  by  Keulegan  and  Carpenter  (1958),  or 

AR  =  C3  co s(30  -  <p3)  -+-  C3  cos(5 0  —  <£5)4 -  . . .  +  Cn  cos (nO  -  <frn),  (1-3) 

where  6  =  ^ .  The  use  of  odd  harmonics  was  justified  since  the  significant  force 
components  occurred  at  the  fundamental  frequency  and  the  3rd  and  5th  harmonics. 
The  force  transfer  model  proposed  by  Keulegan  and  Carpenter  (1958)  and  Sarpkaya 
(1981)  is  then  the  sum  of  Equations  (1.1)  and  (1.3),  or 

/  =  Cdp—  \u\u  -r  CmP—j—a  +  C3  cos(30  -  4>3)  +  C5  cos(50  -  <f>5),  (1.4) 

where  u  ~  -U  cosO,  i.e.  the  flow  is  periodic  with  a  maximum  velocity  equal  to 
f’.  Although  this  model  reduces  the  rms  residual  force  by  approximately  60%, 
it  has  two  limitations-  (1)  the  four  parameters  which  define  the  correction  term 
C3.  C5.  <i> 3  and  <t> 5  have  little  physical  significance,  and  (2)  the  flow  is  assumed  to 
be  periodic.  The  experiments  performed  by  Keulegan  and  Carpenter  (1958)  and 
Sarpkaya  (1981)  prove  that  AR  is  not  the  result  of  the  complications  introduced  in 
wavy  flow.  e.g.  the  limitations  of  the  wave  theory  employed  to  estimate  the  fluid 


particle  kinematics,  the  orbital  fluid  particle  motion,  etc.  Therefore,  attempts  to 
correct  for  this  error  by  adjusting  the  wave  theory  model  may  not  lead  to  significant 
improvement. 

Bird  and  Mockros  (1986)  performed  a  series  of  28  tests  in  which  an  instru¬ 
mented  cylinder  was  accelerated  or  decelerated  in  still  water  in  order  to  estimate 
the  relative  magnitude  of  the  history  force.  Four  cases  were  examined:  (1)  acceler¬ 
ation  from  rest  to  a  constant  velocity,  (2)  deceleration  from  a  constant  velocity  to 
rest,  (3)  acceleration  from  one  constant  velocity  to  a  higher  constant  velocity  and 
(4)  reversal  from  a  constant  velocity  to  a  similar  constant  velocity  in  the  opposite 
direction.  The  results  verify  the  findings  of  Sarpkaya  (1981)  and  suggest  that  the 
relative  magnitude  of  the  history  force  is  significant,  and  at  times  equal  to  the  max¬ 
imum  added  mass  force  during  acceleration  and  over  half  the  steady  state  drag  force 
during  deceleration  for  the  specific  cases  studied.  The  method  employed  by  Bird 
and  Mockros  (1986)  to  estimate  the  history  force  was  somewhat  arbitrary  since  it 
was  presumed  that  the  Morison  equation,  with  constant  coefficients,  can  correctly 
predict  the  instantaneous  drag  and  inertia  components.  Bird  and  Mockros  proposed 
that  the  introduction  of  a  history  force  term  may  account  for  an  important  feature 
of  the  fluid  dynamics  and  reduce  the  variability  between  measured  and  calculated 
forces.  The  force  subdivision  approach  was  based  on  a  discussion  given  by  Hamilton 
(1972).  Hamilton  suggested  that  the  hydrodynamic  force  on  a  submerged  body  in 
nonperiodic  flow  should  be  separated  into  three  parts: 

(1)  the  conventional  drag  component, 

(2)  the  conventional  inertia  component,  and 

(3)  an  explicit  history  term. 
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C.  Regression  Analysis 

Given  the  time  varying  force  acting  on  a  cylinder  and  the  corresponding  fluid 
particle  kinematics,  force  transfer  coefficients  may  be  estimated  in  a  variety  of  ways. 
The  approach  employed  by  most  researchers  in  the  past  has  been  ordinary  least 
squares  regression,  e.g.  Reid  (1958),  Aagaard  and  Dean  (1969),  and  Wheeler  (1970). 
However,  new  statistical  techniques  are  now  available  which  are  particularly  suitable 
for  time  series  data  (Newton,  1988).  Therefore,  a  brief  review  of  regression  analysis 
and  an  alternative  model  fitting  technique  are  presented. 

Regression  analysis  is  the  fitting  of  an  equation  to  a  set  of  values.  The 
equation  predicts  the  response  vector  y  from  a  function  of  the  regressor  matrix 
X  and  parameter  vector  B ,  adjusting  the  parameters  such  that  a  measure  of  fit  is 
optimized.  The  method  used  to  estimate  the  parameters  is  to  minimize  the  sum 
of  the  squares  of  the  differences  between  the  actual  and  predicted  responses.  In 
matrix  notation 

y  —  XB  +  e,  (1.5) 

such  that  B  minimizes 

S(B)  =  (y  -  XB)r(y  -  XB).  (1.6) 

This  analysis  is  based  on  several  assumptions,  including:  (1)  that  the  expected  value 
of  the  errors  (e)  is  zero  and  (2)  that  the  errors  are  uncorrelated  across  observations. 
When  regression  is  performed  on  time  series  data,  the  errors  are  often  autocorre- 
lated,  violating  the  second  assumption  and  perhaps  the  first.  In  this  case  there 
are  several  new  methods  to  estimate  B.  One  such  method  is  regression  with  auto¬ 
regressive  (AR)  errors.  This  is  an  iterative  procedure  in  which  the  following  steps 
are  performed: 


(1)  use  ordinary  least  squares  to  find  initial  estimates  of  B  and  e, 

(2)  fit  an  AR  model  to  e, 

(3)  create  a  new  response  vector  z  and  regressor  matrix  W  by  applying  the  AR 
filter  of  step  (2)  and 

(4)  return  to  step  (1)  using  2  and  W  instead  of  y  and  X . 

This  process  is  repeated  until  successive  iterations  result  in  the  same  value  of  AR 
order.  A  more  thorough  discussion  along  with  the  software  required  to  perform  the 
above  analysis  is  given  by  Newton  (1988). 

D.  U.S.  Naval  Academy  Experimental  Program 

In  order  to  gain  insight  into  the  applications  and  limitations  of  the  relative 
motion  form  of  the  Morison  equation,  a  multiple  phase  experimental  program  was 
conducted  at  the  U.S.  Naval  Academy,  Hydromechanics  Laboratory  (Shields  and 
Hudspeth,  1985).  The  results  from  the  earlier  phases  are  reported  by  Dawson 
(1984)  and  Dawson,  Wallendorf  and  Hill  (1986).  The  data  analyzed  in  this  thesis 
is  the  result  of  an  extension  of  these  experiments  to  include  both  relative  motion 
and  random  wave  loading.  The  USNA  test  assembly  was  subjected  to  random 
waves  of  approximate  Bretschneider  spectra  with  a  significant  wave  height  of  20 
cm  and  a  dominant  wave  period  of  2.38  seconds.  The  maximum  Reynolds  number 
was  approximately  2  x  104  and  the  maximum  Keulegan- Carpenter  number  was 
approximately  50. 

The  test  assembly  consisted  of  a  smooth  stainless  steel  rod  of  one  meter  length 
and  2.54  cm  diameter.  The  cylinder  is  attached  to  a  free-to-surge  aluminum 
subcarriage  which  is  supported  on  two  linear  bearings  and  restrained  be  elastic 
springs.  The  test  assembly  is  illustrated  in  Figure  1.  Instrumentation  used  in  the 
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experiments  included  a  high  frequency  sonic  transducer,  a  variable  reluctance  force 
gage  and  a  resistance  wave  staff.  Instrument  signals  were  passed  through  an  analog 
filter  with  20  Hz  cut-off  frequency  and  an  analog-to-digital  converter  before  storing 
on  computer  disk.  The  acquisition  rate  was  51.2  samples  per  second.  The  accuracy 
of  the  measured  surface  elevation,  surge  displacement,  and  force  are  estimated  to 
be  0.02  cm,  0.02  cm,  and  0.05  N,  respectively  (Dawson,  Wallendorf  and  Hill,  1986). 
The  data  is  analyzed  in  the  next  chapter  using  time  series  analysis  techniques. 
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CHAPTER  II 
DATA  ANALYSIS 

A.  Fluid  Particle  Kinematics 

The  first  step  in  analyzing  the  USNA  data  is  to  estimate  the  fluid  particle 
kinematics  from  the  measured  surface  elevation  time  series.  This  is  accomplished 
using  harmonic  analysis.  The  measured  surface  elevation  is  represented  as  an  infinite 
sum  of  small  amplitude  waves,  with  closely  spaced  frequencies,  and  random  phase 
angles 

OO 

rj(x,t)  =  ^  Ancos(knx  -  u>nt  +  <£n),  (2.1) 

n—  1 

where  77  denotes  the  elevation  of  the  water  surface  above  the  mean  water  level,  An  is 
the  amplitude  of  the  nth  wave  component,  and  fcn,  ivn,  and  <f>n  are  the  corresponding 
wave  number,  frequency,  and  phase  (Borgman,  1972).  The  amplitudes  An  and  phase 
angles  4>n  are  calculated  most  efficiently  using  the  fast  Fourier  transform  (FFT) 
algorithm.  The  horizontal  velocity  and  acceleration  at  (x,z,t)  are  then  given  by 

.  V"'  coshfcn(h  +  z)  .  _ 

u(x,z,t)  =  - sinh  h  '4n  cos^x  -a,"*  +  ^n),  (2.2) 

n=  1  n 

and 

,  ,  ,  cosh  fcn(h  +  z)  ,  .  . ,  ,  .  . _ _ 

a{x,z,t)  =  2]  - sinhy  h — -  An  sm(A;x  -  wnt  +  <f>n),  (2.3) 

n=  1  n 

where  the  wave  number  and  frequency  are  related  through  the  dispersion  relation 

=  tanh  knh.  (2-4) 

In  Equations  (2.1)  through  (2.3)  the  Fourier  components  are  summed  from 
n  =  1  to  00,  however,  wave  records  generally  consist  of  discrete  samples.  A  sampling 


interval  of  1/1  Oth  to  1  /20th  of  the  dominant  wave  period  is  often  recommended. 
The  highest  frequency  that  can  be  used  in  the  Fourier  summation  is  then 


(2.5) 


where  At  is  the  sampling  interval.  This  frequency  is  known  as  the  Nyquist. 
frequency.  The  dominant  wave  period  of  the  USNA  data  is  2.38  seconds.  A 
sampling  interval  of  approximately  1/1 5th  of  the  dominant  wave  period  was  chosen, 
or  At  =  0.15624  seconds.  Therefore,  the  Nyquist  frequency  is  approximately  3.2 
Hz.  This  sampling  interval  reduces  the  amount  of  data  to  be  processed  without 
loss  of  any  significant  high  frequency  components.  The  spectral  density  estimates 
before  and  after  sampling  were  compared  to  ensure  that  the  measured  spectra  were 
represented  accurately. 


B.  Regression  Model  Formulation  and  Solution 

Given  the  time  varying  force  acting  on  the  cylinder,  and  the  corresponding 
fluid  particle  kinematics  and  cylinder  dispacement,  force  transfer  coefficients  can 
be  estimated  using  the  regression  with  AR-errors  technique  described  in  Chapter  1. 
The  total  force  acting  on  the  cylinder  at  time  t  may  be  considered  as  the  linear 
combination  of  three  terms:  (1)  the  drag  component  F d  given  by 

D 

Fd(1)  =  p-y  J  h \u{z,t)  -  x(t)\{u{z,t)  -  x(t)}dz,  (2.6) 

(2)  the  inertia  component  Fj  given  by 

Fi{i)  —  p  r  /  {a(z,<)  -  x{t)}dz,  (2.7) 

J  —  h. 

and  (3)  a  force  component  resulting  from  the  acceleration  of  the  fluid  by  the  cylinder, 
Ff  given  by 

FF(t)  =  pVi(t),  (2.8) 


where  h  is  the  submerged  length  of  the  cylinder  and  V  is  the  volume  of  displaced 
fluid.  The  last  term  may  be  combined  with  the  measured  force  acting  on  the  cylinder 
and  denoted  by  F(t)  to  give 


F(t)  =  CdFD(t)  +  CmFj(t).  (2.9) 

Equation  (2.9)  is  now  in  the  form  of  the  regression  model  y  —  XB  +  e,  where 
X  =  ( Fd,Fj )  and  B  =  {Cd,C7n)T .  The  best  linear  unbiased  estimate  of  B  results 
in: 

Cd  =  0.93  and  Cm  =  1.73. 

A  thorough  time  series  analysis  of  the  errors  results  in  the  following  observations: 
(1)  the  errors  e  are  highly  autocorrelated,  (2)  the  cross-correlations  between  the 
prewhitened  error  time  series  and  the  prewhitened  cylinder  velocity  and  surface 
elevation  time  series  are  relatively  small  (less  than  0.4),  and  (3)  the  cross-correlation 
between  the  prewhitened  error  time  series  and  the  prewhitened  measured  force  F  is 
large  (approximately  0.8  at  lag  zero).  Based  on  these  observations  one  may  conclude 
that  the  errors  are  not  the  result  of  the  extension  of  the  Morison  equation  to  include 
relative  motion,  but  are  indicative  of  a  systematic  modeling  error  in  the  original 
Morison  equation.  The  primary  source  of  this  error  is  determined  in  the  following 
section. 

C.  The  Best-Fit  Relative  Motion  Morison  Equation 

Figure  2  compares  the  best-fit  relative  motion  Morison  equation  with  the 
measured  force  over  a  16  second  realization.  The  Morison  equation  “tracks” 
the  measured  time  series  fairly  well,  but  underpredicts  the  rms  force  by  21%. 


y.\ 


Relative  Motion  Morison  Equation.  Comparison  of  Measured  and 
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The  multiple  correlation  coefficient,  defined  as 


R 2  =  1  - 


(2.10) 


where  and  cr^  are  the  sample  variance  of  the  errors  and  the  measured  force 
realizations,  is  0.77.  One  can  see  qualitatively,  that  the  predicted  time  series  is 
“smooth”,  while  the  measured  time  series  is  somewhat  “wiggly”.  This  is  due  to  the 
presence  of  a  series  of  high  frequency  force  components  that  the  Morison  equation 
lacks.  This  discrepancy  is  illustrated  using  the  standardized  spectral  density,  defined 


Sn(o>) 


k  1 


versus  u>  = 


(2.11) 


for  k  =  l,...,[y]  +  l,  where  Su  is  the  usual  two-sided  spectral  density,  a 2  is  the 
sample  variance,  and  n  is  the  total  number  of  observations  (Newton,  1988).  The 
standardized  two-sided  spectral  density  estimates  for  both  predicted  and  measured 
forces  are  shown  in  Figure  3. 

The  high  frequency  residual  force  is  readily  apparent,  with  a  peak  at  approxi¬ 
mately  five  times  the  fundamental  frequency.  The  cause  of  this  error  can  be  deter¬ 
mined  through  cross-spectra  analysis. 

Recall  that  the  drag  and  inertia  components  serve  as  the  independent  vectors 
in  the  regression  model,  while  the  measured  force  is  the  dependent  vector.  Figures 
4  and  5  show  the  cross-spectra  of  the  independent  and  dependent  vectors. 

Figure  4  is  a  graph  of  the  coherency  spectrum  Wiz(u))  which  measures  the  rela¬ 
tionship  between  the  amplitudes  of  the  sinusoids  in  the  two  univariate  realizations 


at  frequency  u>.  The  coherence  is  given  by 


W12{U)  = 


\SnH\ 


(2.12) 
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Coherency  Spectra:  Independent  &  Dependent  V 
Occur  at  Both  the  Fundamental  and  5th  Harmonic 
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where  Sn  and  522  are  the  usual  two-sided  univariate  spectral  densities,  and  S n  is 
the  complex  cross-spectral  density.  Figure  5  is  a  graph  of  the  phase  spectrum  <£12(u>) 
which  measures  how  out  of  phase  the  frequency  components  for  the  univariate 
realizations  tend  to  be.  The  phase  is  given  by 

^2(a,)  =  tan-1-{^^},  (2.13) 

c12(u;) 

where 

Ci2(w)  =  Re  (512(w)}  and  q12(uj)  =  -Im{5i2(u;)}  ,  (2.14) 

are  the  cospectral  density  and  the  quadrature  spectral  density,  respectively. 

Two  important  observations  can  be  made  from  the  cross-spectra:  (1)  the 
coherence  is  significant  at  both  the  fundamental  and  5th  harmonic,  and  (2)  a 
considerable  phase  shift  exists  between  the  fundamental  and  the  5th  harmonic. 
Since  the  relative  motion  Morison  equation,  which  contains  only  instantaneous  and 
time  averaged  parameters,  is  able  to  capture  the  fundamental  frequency  component, 
but  not  the  5th  harmonic,  one  may  conclude  that  the  5th  harmonic  lags  the 
fundamental  frequency  component.  That  is,  the  5th  harmonic  is  the  result  of  flow 
history.  If  this  is  true,  then  considerable  improvement  in  force  prediction  can  be 
achieved  through  the  use  of  a  force  transfer  model  containing  lagged  versions  of 
the  independent  vectors,  i.e.  history  terms.  Such  a  model  is  presented  in  the  next 
section. 

D.  A  Force  Transfer  Equation  with  Memory 

Recall  that  the  regression  model  used  to  estimate  the  best-fit  relative  motion 
Morison  equation  is 

F(t)  =  CdFD(t)  +  CmFr{t).  (2.15) 
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Instead  of  the  usual  approach  of  using  time  averaged  force  transfer  coefficients  Cd 
and  Cm,  consider  F(t )  as  the  sum  of  a  filtered  version  of  the  drag  component  time 
series  plus  a  filtered  version  of  the  inertia  component  time  series,  or 

oo  oo 

F(t)  =  £  <*kFD(t-k)+  ]T  faF^t-k),  (2.16 a) 

k=  —  o o  fc=  — oo 


with  the  following  constraints 


f  0,  if  k  <  0; 
|  0,  if  k  >  l 


and  /3k  = 


0,  if  k  <  0; 
0,  if  k  >  l 


(2.166) 


where  a  and  f3  are  force  transfer  coefficient  vectors,  and  /  is  the  maximum  lag  used 
in  the  model,  i.e.  I  is  a  measure  of  the  memory.  Note  that  when  Z  =  0  Equation 
(2.16)  reduces  to  Equation  (2.15),  the  relative  motion  Morison  equation.  Thus,  the 
instantaneous  force  acting  on  the  cylinder  is  now  represented  by  the  superposition 
of  a  filtered  version  of  the  drag  component  and  a  filtered  version  of  the  inertia 
component,  both  containing  only  non-negative  lags.  In  this  thesis,  Equation  (2.16) 
will  be  called  the  “force  transfer  filters”.  The  improvement  in  force  prediction,  and 
the  tradeoff  between  model  complexity  and  goodness-of-fit  is  illustrated  in  the  next 
section. 


E.  Model  Complexity  vs.  Goodness-of-Fit 

Figures  6  through  10  illustrate  the  improvement  in  force  prediction  with 
increasing  memory.  As  the  memory  is  increased,  the  error  between  rms  forces  is 
reduced  while  the  multiple  correlation  coefficient  R 2  is  increased.  These  calcualtions 
were  performed  on  27  realizations  containing  a  total  of  432  seconds  of  data 
(approximately  180  waves).  As  is  the  case  with  most  regression  models,  there 
is  a  tradeoff  between  model  complexity  and  goodness-of-fit.  For  the  USNA  data  a 
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ten  lag  memory  was  chosen  as  the  optimal  model.  Using  the  force  transfer  filters 
with  /  =  10  reduces  the  mean  rms  error  from  24%  to  9%  while  increasing  the  mean 
multiple  correlation  coefficient  from  0.71  to  0.83.  Since  the  sampling  interval  is 
At  =  0.15624  seconds,  the  10  lag  memory  corresponds  to  1.5624  seconds.  This  can 
be  nondimensionalized  by  dividing  by  the  dominant  wave  period  T  =  2.38  seconds 
to  give 

l 

At-  =  0.66. 

T 

Thus,  the  optimal  model  for  the  USNA  data  contains  a  memory  equal  to  two  thirds 
of  the  dominant  wave  period. 

In  addition  to  decreasing  the  rms  residual  force  by  an  average  of  170%,  the 
force  transfer  filters  capture  the  5th  harmonic  not  accounted  for  in  the  Morison 
equation.  These  small  amplitude,  high  frequency  force  components  may  govern  the 
fatigue  analysis  of  offshore  structures.  As  shown  in  Figure  11,  the  ten  lag  memory 
force  transfer  filters  duplicate  the  measured  spectral  density  almost  exactly. 


;ag  Memory.  (Morison  equation.)  Comparison  of  Measured 
Realizations  -  RMS  Error=21%,  R2=  0.77. 


y.  Comparison  of  Measured  and  Predicted  Force  Realizations 
Transfer  Filters  -  RMS  Error=13%,  R2  —  0.82. 


Memory.  Comparison  of  Measured  and  Predicted  Force  Realizations 
;  Force  Transfer  Filters  -  RMS  Error=6%,  f?J=0.86. 


26 


i 

ft'*-' 

* 


✓ 

> 


% 


CHAPTER  III 
NUMERICAL  SIMULATION 

A.  Tension  Leg  Platform  Model 

The  purpose  of  this  chapter  is  to  provide  a  simple  example  of  the  application  of 
the  force  transfer  models  discussed  in  Chapter  2.  A  tension  leg  platform  is  modeled 
as  a  nonlinear,  single  degree  of  freedom  system  (SDOFS)  by  considering  motion  only 
in  the  surge  direction.  Riser  and  tendon  dynamics  are  neglected.  Furthermore,  it 
is  assumed  that  diffraction  and  reflection  of  the  waves  are  negligible.  The  equation 
of  motion  is  established  by  equating  the  sum  of  the  inertia,  damping,  stiffness  and 
external  forces  at  time  t*: 

mii  +  cii  +  kx  =  Fi,  (3-1) 

where  m  is  the  mass  of  the  platform  and  deck  equipment,  c  is  the  damping,  k  is 
the  stiffness  and  Fi  is  the  hydrodynamic  force  resulting  from  the  fluid-structure 
interaction. 

Following  the  derivation  by  Malaeb  (1982),  the  stiffness  is  obtained  by  displac¬ 
ing  the  platform  an  arbitrary  distance  in  the  surge  direction  (Figure  12). 

Summing  the  horizontal  forces  gives 

kx  =  N(To  +  AT1)  sin  if)  (3-2) 

where  k  is  the  stiffness  in  the  x-direction  due  to  an  arbitrary  surge  displacement,  N 
is  the  number  of  tension  legs,  T0  is  the  pretension  of  each  leg  and  A T  is  the  increase 
in  tension  per  leg,  given  by 


i 


AT=kcAd, 


(3-3) 
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where  kc  is  the  equivalent  stiffness  per  leg.  In  most  cases  AT  is  less  than  20%  of 
the  pretension  T0.  The  change  in  cable  length  Ad  and  the  angle  of  inclination  tp 


are  given  by 


Ad  =  yd2  4-  22  —  d, 


sin  ip  ~ 


>d?  + 


Substituting  equations  (3.3)  through  (3.5)  into  Equation  (3.2)  gives 


k  =  (To  +  kcy/  d?  +  x 2  —  kcd)- 


U2  +  x2 


Thus,  the  restoring  force  resulting  from  an  arbitrary  surge  displacement  is  a  function 
of  the  magnitude  of  the  displacement  as  well  as  the  material  properties,  pretension, 
and  elongation  of  the  tendons.  The  form  of  Equation  (3.6)  is  that  of  a  “hardening 
spring”.  When  the  surge  displacement  is  small  (*  <<  d )  Equation  (3.6)  can  be 
reduced  to  a  linear  function  of  x: 


The  restoring  force  kx  vs.  x  is  plotted  in  Figure  13  using  both  Equations  (3.6)  and 
(3.7)  with  kc  —  7.0(107)  N/m,  To  =  4.8(107)  N,  and  d  =  200  m.  The  data  listed  in 
Table  1  are  used  for  simulation  purposes. 
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Table  1.  Typical  TLP  Data  Based  on  K:rk  and  F.tok  (1979) 


ft 

— 

Structural  Component 

Dimensions 

M; 

Diameter  of  4  corner  columns 

]  6  in 

Diameter  of  4  middle  columns 

1 

3.5  m 

“V** 

Diameter  of  2  cross  braces 

6  m  j 

■v* 

Spacing  of  corner  columns 

70  in 

Draft,  h 

35  m 

,v 

Total  mass  in  air.  m 

5.4(10' )  kg 

Fluid  added  mass  in  surge 

3.3(10")  kg 

i 

Cable  stiffness  per  leg,  kc 

7.0(10")  N  m 

Pretension  of  each  leg.  T0 

107  -  108  N 

,V 

Water  depth 

200  -  1000  m 

>/ 

i 

Critical  damping  ratio,  £ 

0.00  -  0.15 

B.  Solution  Procedure 
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The  nonlinear  equation  of  motion  can  be  solved  using  one  of  several  direct 
integration  procedures.  The  Newmark-i?  method  is  chosen  because  of  its  well  known 
convergence  and  stability  criteria  (Newmark,  1962).  The  Newmark-i  equations  are 
used  to  determine  the  velocity  and  displacement  of  the  structure  at  time  t,-_j  based 
on  the  corresponding  value-  at  time  and  the  acceleration  at  time  <,-i-  Since  the 

The  advantage  of  using  such  a  procedure  is  that  any  nonlinearities,  such  a.-  the 

I  LP  stiffness  given  by  Equation  (3.6),  can  easily  be  included  in  the  analysis.  The 
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Newmark-/?  equations  are  as  follows: 


*t+i  =  i.  +  (1  -  7)*»A<  +  7*i+i  At 


xi+i  =xi  +  ijAf  +  (0.5  -  /3)xi{At)2  +  pxi+1(At)2  (3.9) 

where  7  is  taken  as  0.5  to  avoid  spurious  numerical  damping  and  f3  is  taken  as  0.126 
to  ensure  numerical  stability  with  convergence  (Newmark,  1962).  The  assumed 
value  of  ij4-i  is  compared  with  the  corresponding  value  calculated  from  the  equation 


of  motion: 


*t+l  —  (-fi+l  c*»+ 1 

m 


(3.10) 


Equations  (3.8)  through  (3.10)  are  solved  in  an  iterative  manner  until  the  assumed 
and  calculated  values  of  ij+1  converge.  The  process  is  then  repeated  for  subsequent 
time  steps. 

In  order  to  ensure  numerical  stability,  a  two-pass  procedure  is  used  wherein 
the  response  of  the  TLP  is  first  estimated  using  the  relative  motion  form  of  the 
Morison  equation.  The  drag  and  inertia  force  components  are  then  filtered  and 
superimposed  using  Equation  (2.16)  to  obtain  an  adjusted  force  time  series.  This 
time  series  is  then  used  as  input  to  the  Newmark-/?  algorithm,  thus  obtaining  an 
adjusted  response  time  series  (Figure  11). 


C.  Simulation  Results 


The  results  from  two  typical  numerical  simulation  cases  are  presented  in  Figures 
15  through  24.  In  case  one,  the  random  wave  profile  is  scaled  such  that  the  maximum 
wave  height  is  approximately  10  meters.  In  case  two,  the  same  wave  profile  is  used, 
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but  scaled  such  that  the  maximum  wave  height  is  30  meters.  The  force,  surge, 
and  tension  time  series  are  plotted  after  both  the  first  and  second  pass  of  the 
Newmark-/?  algorithm.  In  both  cases,  the  force  transfer  filters  (2nd  pass)  result  in 
a  substantial  increase  in  high  frequency  force  components.  These  components  are 
more  pronounced  when  the  wave  height  is  small  and  the  wave  period  is  short.  This 
may  be  due  to  the  phase  shift  resulting  from  the  spacing  of  the  columns.  If  the  vave 
length  is  equal  to  twice  the  column  spacing  (140  meters),  then  the  resulting  forces 
on  opposing  columns  will  partially  cancel  with  one  another.  In  all  cases,  the  surge 
response  is  smooth  because  of  the  large  inertia  of  the  TLP.  Since  the  variation  in 
tension  was  represented  as  a  function  of  the  surge  displacement  by  neglecting  the 
dynamics  of  the  tension  legs,  the  tension  time  series  is  also  smooth.  If  the  dynamics 
of  the  tension  legs  were  included  by  modeling  the  TLP  as  a  MDOFS,  then  these 
high  frequency  force  components  should  also  occur  in  the  tension  time  series.  Such 
an  analysis  could  explain  the  tension  leg  “ringing”  phenomenon  observed  in  the 
Hutton  TLP  model  tests  (Mercier,  1982). 

Although  both  cases  produce  similar  results,  the  effect  of  the  force  transfer 
filters  is  more  pronounced  in  case  two  where  the  maximum  wave  height  is  much 
greater  than  in  case  one.  The  principle  limitation  of  the  use  of  these  filters  is 
that  the  scale  effect  is  unknown  at  this  time.  Like  the  Morison  equation,  the 
force  transfer  filters  should  be  calibrated  for  a  specific  design  application,  taking 
into  consideration  the  wave  theory  employed  to  estimate  fluid  particle  kinematics, 
the  ocean  wave  design  spectra,  and  the  size,  shape,  and  interaction  of  the  various 
structural  components  which  comprise  the  platform. 
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CHAPTER  IV 


CONCLUSIONS 


In  this  study,  experimental  data  were  analyzed  to  gain  insight  into  the  appli¬ 
cations  and  limitations  of  the  relative  motion  form  of  the  Morison  equation  for  the 
prediction  of  hydrodynamic  forces  on  a  free-to-surge  vertical  cylinder.  A  new  force 
prediction  procedure  was  presented  to  account  for  the  effects  of  flow  history  on  the 
instantaneous  force.  The  following  conclusions  are  made  based  upon  this  research: 

(1)  Small  scale  model  tests  suggest  that  the  effect  of  flow  history  is  significant  in 
the  prediction  of  hydrodynamic  forces. 

(2)  The  relative  motion  form  of  the  Morison  equation  can  not  duplicate  the  high 
frequency  force  components  measured  in  small  scale  model  tests, 

(3)  In  this  study,  the  inclusion  of  explicit  history  terms  captured  the  high  frequency 
force  components,  increased  the  multiple  correlation  coefficient  from  0.71  to 
0.83,  and  decreased  the  root-mean  square  error  from  24%  to  9%. 

(4)  The  procedure  to  include  history  terms  can  be  incorporated  into  the  response 
simulation  of  more  complex  offshore  structures  such  as  the  TLP. 

The  results  from  this  research  suggest  that  significant  improvements  in  hy¬ 
drodynamic  force  prediction  can  be  obtained  by  including  explicit  history  terms. 
A  strong  research  effort  is  required  to  determine  the  scale  effect  on  these  history 
terms.  Future  research  should  include  both  large  scale  model  tests  and  field  data 
analysis.  The  results  may  prove  to  be  of  great  importance  in  the  fatigue  analysis  of 
offshore  structures  and  may  explain  some  of  the  variability  between  measured  and 
calculated  forces. 
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APPENDIX  A 

NOMENCLATURE 

amplitude  of  nth  wave  component 
horizontal  fluid  particle  acceleration 
regression  parameter  vector 
added  mass  coefficient 
drag  coefficient 
inertia  coefficient 

remainder  force  Fourier  coefficients 
structural  damping  coefficient 
cospectral  density 
diameter 

length  of  tension  leg  at  zero  surge  displacement 

change  in  length  of  tension  leg 

total  hydrodynamic  force 

hydrodynamic  drag  component 

hydrodynamic  inertia  component 

force  resulting  from  the  acceleration  of  fluid 

hydrodynamic  force  pe"  unit  length 

Nyquist  frequency 

acceleration  of  gravity 

water  depth 

structural  stiffness 

equivalent  stiffness  per  tension  leg 
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X 

X ,  X,  X 

y 


wave  number  of  nth  wave  component 

maximum  lag  used  in  a  filter  i.e.  memory  length 

mass  of  structure 

number  of  tension  legs 

quadrature  spectral  density 

multiple  correlation  coefficient 


*■ 

AR 

remainder  force  i.e.  error 

r- 

T 

wave  period 

To 

pretension  of  each  tension  leg 

& 

AT 

increase  in  tension  per  tension  1< 

-■  , 

t 

time 

!j 

At 

change  in  time 

S(B ) 

sum  of  the  square  of  the  errors 

a 

■Siii  £22 

univariate  spectral  density 

s 

c 

complex  cross-spectral  density 

, 

u 

horizontal  fluid  particle  velocity 

Hri2 

coherency  spectrum 

regressor  matrix 

surge  displacement,  velocity,  and  acceleration 

response  vector  used  in  regression  analysis 

elevation  from  mean  water  level 

drag  filter  coefficients 

Newmark-/?  parameter 

inertia  filter  coefficient 

Newmark-/?  parameter 
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error  time  series 

surface  elevation  time  series 

nondimensionalized  time  defined  as  ^ 

3.141 ... 

fluid  mass  density 
sample  variance 

phase  angle  of  the  nth  wave  component 
phase  spectrum 

angle  of  inclination  in  the  xz  plane 
frequency  given  by  for  fc  =  l,...,[j]  +  l 
frequency  of  the  nth  wave  component 
volume  of  displaced  fluid 
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APPENDIX  B 

TIME  SERIES  ANALYSIS  SOFTWARE 

The  statistical  analysis  results  reported  in  this  thesis  were  obtained  using 
TIMESLAB.  This  computer  program  is  essentially  a  time  series  analysis  language 
consisting  of  approximately  150  commands.  The  software  is  provided  with  the 
text:  TIMESLAB:  A  Time  Series  Analysis  Laboratory  (Newton,  1988).  Some  of 
the  features  which  were  used  in  this  research  include:  regression  with  AR-errors, 
a  test  for  white  noise,  auto-  and  cross-correlation,  coherence  and  phase  spectrum 
estimation,  nonparametric  spectral  density  estimation,  and  a  test  for  statistical 
independence  of  two  univariate  realizations. 

The  TIMESLAB  user  can  create  files  or  “macros”  containing  commands. 
The  macros  used  to  estimate  and  plot  the  best-fit  relative  motion  form  of  the 
Morison  equation  and  the  best-fit  force  transfer  filters  are  provided  below.  Note 
that  FILTERS. MAC  calls  another  macro,  REGAR.MAC,  which  is  not  listed  here. 
REG AR. MAC  is  provided  as  part  of  the  TIMESLAB  software  package. 
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FILTERS. MAC 


I nput  : 


PAUaE 
: start 


TIMESLAE  macro  to  estimate  and  plot  tne  oest-f't  reiat-ve 
motion  Morison  equation  and  tie  best-fit  force  transfer 
filters  using  regression  with  AR-errors. 


yi 

X  1 

x2 

max ' ag 


total  force 
drag  component 
inertia  component 
time  w i nacw 
max i mum  1 ag 


n_.  ;  regress’On  parameters 

ma xp  =  30 

; z  ,  ma (c  'e  m i . s« i P  1 
maxc  =  .Ti  ■ 

;  Sk  l p 


y  =  EXTRACT(y1  .  l  . j  ) 
x11ag=EXTRACT(x1 ,  i  ,  j  ) 
x2 1 ag  =  EXTRACT ( x2 .  i  .  j  ) 
X  =  <x 1 1 ag , x2 1 ag> 


instantaneous  total  force 
instantaneous  drag  force 
instantaneous  'inertia  force 
form  regressor  matrix 


m  =  2 

MACRO ( rega- . start  ) 

_ ABEL ( yy )=' Mor ■ son  Eduat i on 
va  xMI M(  -./max',  i max , , r i r i .  iminl 

MAXM I N ( y>  ,  n .  ymax2  ,  1  max  .  yrr,  l  n2  .  im  i  n ) 
yaes  =  ABS ( ym i " *  1 
ysave=y 
y=y+yabs+ymax 1 
vsave  i  =y 
_  A  3  E  L  (  y  )  - 

MAXMINt  y  ,  n  .  /max'  ,  •  r  a  *  .  .-rim  .  irom) 
Pl0TS:ZE(480.  1 20 . 55 , 3C . S ,  10,8.  1,0,01 
t ime* LINE ( 101 . -0. 1562 A .0- 15624) 

MACRO (error, start) 

LABEL ( t me ) x ' T ime  (sec).  Error=#errf%; 
PL0T0N 

PLOT ( 1 1  me , y , n , 0. 1 6 , ym l n2 . ymax 1 )  ; 

;  L0T( t  ime , yy . n ,0, 16 , ymin2 , ymax 1 )  ; 

plotoff 

Cd=beta [ 1 ] 

Cm=beta[ 2 ] 

LIST! Cd , Cm )  : 

3R0MPT0FF 

m  =  rn  1 
k  =  ' 

1  =  J 

km i n=k -max ’ ag 

;  1  OOP  • 


•  regress  on  instantaneous  values 


R2=^R2P' 

plot  measured  and  predicted  force 
using  the  Mori  son  equation 

list  drag  and  inertia  coefficients 


loop  to  get  lagged  values 


k  =  k  -  1 
1  =  1-1 
' ag- i -k 

’  C  -  {  '  3  2  I 

*  '  "  • 
x  2  1  ag  -  E.‘  ”PACV  i  x  2  ,  k  .  i  ) 
; =<*  lag,x2 :ag> 

IF  (k  .eq.kfnm,  stop  ) 

GOTO  (loop) 

:  Stop 


form  '•egressor  matrix 
using  lagged  values 


end  of  loop 
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CLEANlxilag.x21ag,k,l ) 
MACRO! "egar .start ) 


regress  with  history  te-ns 


_  AE  ;  l  I  yy  !  -  ■  ►  tf  with  xlagx  Lao  Meniere-' 

,  =  .save' 

-AEE. (  >  i= 

MA  >  M  r  N (  yy . n . ymaxC ,  i max . ym 1 n2 . 1  mi n I 
MACRO! error . start ) 

-A£E_(  1 1  me  )  =  '  T  i  me  (sec).  Error  =  «err«}:;  R2=?-R2ca 
p-_0t0N 

~  “  1  t  i  me  .  y  .  n ,  w  .  16  .  ym  i  n2  .  yrr.ax  i  1  ;  p'ct  measure  anc  DreO'cted  ;orte 

*  *  i  me  ,  yy  ,  r. .  C .  1 6  ,  ym  l  n2  .  ymax  1  )  :  using  the  force  transfe'  fi’te-' 

P-CTO-f 

p.c'size; c ) 


ERROR  Macro  to  compute  the  rms  error  and  multiple  correlation 
coef f i c i ent 


Input:  ysave  =  Measured  Fo^ce 

'/y  =  Predicted 

n  =  length  of  arrays 


PAuSE 
: start 

yrms=MMULT ( ysave . y Save , 1 , n. 1 > 
yr  ms* {yrms/n} *0 . 5 

yyrms  =  MMULT ( yy . yy . 1 . n  ,  1  ) 

, yrrs5 v/yrms  n> *C  .  5 

mrms=yrms [ 1 ] 
prms=yyrms[ 1 ] 
c 1  ear! yrms , yyrms ) 
error* {mrms-prms}/m rms* 100 
err  =  C 

err : =error 

del ta=error-err 

1 F ( del ta . 1 t . 0 . 5 , sk ip ) 

err=err*  1 

:  sk  i  p 

rhc=CORR( res.n, 30.0. 1 . se2 . per ) 
rho* CORR ( y , n . 30 , 0 .  1 . Sy2 . per  ) 

R2  = "  0- see ■ s/2 
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